version 14.2
clear all
set more off
set matsize 11000
set maxvar 10000
set more off

global dir "C:\Users\antoi\Dropbox\Peer and religious schools\Peer Religion Replication"


cd "$dir\Raw_Data"

use base_data.dta, clear
drop _merge
merge m:1 school1 using "Schools12.dta"

/*Drop all the students that started school as Freshman and students who are in
Medical or graduate school*/
keep if FR_nonFR_status =="never enrolled as FR" | FR_nonFR_status ==""
drop if MAJR_CODE=="DMED"
drop if Class=="GR"

drop if Class =="UP"
drop if Class=="SP"
drop if Class=="PG"
drop if Class=="UP"


destring Grade, generate (grade) force
sort Term Subject CourseNum, stable
by Term Subject CourseNum: egen mean_coursegrade=mean(grade) 
by Term Subject CourseNum: egen sd_coursegrade=sd(grade) 
by Term Subject CourseNum: gen st_coursegrade=(grade-mean_coursegrade)/sd_coursegrade

tostring Term, generate(str_geocode)
gen Year = substr(str_geocode,1,4)
destring Year , generate(year) force

sort ID Term, stable
by ID Term: egen SemGPA = mean(st_coursegrade)
sort ID year, stable
by ID year: egen YearGPA = mean(st_coursegrade)

drop CourseNum Subject SECTION_CREDIT_HOURS Grade instructor_sex sd_coursegrade st_coursegrade mean_coursegrade grade dummy_withdrew
duplicates drop 

*First major
sort ID, stable
by ID: stripolate  first_major2 year, gen(firstrealmajor2) groupwise
by ID: stripolate  first_major year, gen(firstrealmajor) groupwise /*Including Majorless*/

*This is for first year major chosen (i.e. not majorless).
gen major_science=1 if firstrealmajor2=="MATA" | firstrealmajor2=="MATS" | firstrealmajor2=="APPM" | firstrealmajor2=="PHYS" | firstrealmajor2=="CMPS" | firstrealmajor2=="CHEM" | firstrealmajor2=="GEOL" 
replace major_science=1 if firstrealmajor2=="CONS" | firstrealmajor2=="EECE" | firstrealmajor2=="EELE" | firstrealmajor2=="ECIV" | firstrealmajor2=="ECMP" | firstrealmajor2=="EMCH"  | firstrealmajor2=="CHME" | firstrealmajor2=="INDE"
replace major_science=0 if firstrealmajor2!="" & major_science==. 

/*Variale bla created below determines the number of semester a student has been in school
at the observation.*/
sort ID Term, stable
by ID : gen bla= _n 
by ID: egen maxbla = max(bla)
by ID: egen minbla = min(bla) if year>=2004
gen newbla = 1 if minbla==bla

*This is to know which was last semester for students
gen lastsem = 1 if bla==maxbla 


*****************************************************************************
*																			*
*				Treatment variables and Controls							*
*																			*
*****************************************************************************
*Create Treatment Variables
*Create year ID for advisor
gen advisorIDyear=fake_advisor_id+year

*Fill in school type for some of the missing years*
sort ID, stable
by ID: egen mmuslim=mean(muslim)
by ID: egen mchristian=mean(christian)

drop if mmuslim==.5 | (mmuslim>1 & mmuslim<98)
drop if mchristian==.5 | (mchristian>1 & mchristian<98)


drop if mmuslim==98 | mchristian==98

gen islam_school=(mmuslim==1)
gen christ_school=(mchristian==1)
gen secular_school=1 if islam_school!= 1 & christ_school!=1
replace secular_school=0 if secular_school==.

sort fake_advisor_id year, stable
by fake_advisor_id year: gen newtotal=_N -1

by fake_advisor_id year: egen sum_islam=sum(islam_school)
by fake_advisor_id year: gen summus1out=sum_islam-islam_school
by fake_advisor_id year: gen prop_islam=summus1out/newtotal

by fake_advisor_id year: egen sum_christ=sum(christ_school)
by fake_advisor_id year: gen sumchris1out=sum_chris-christ_school
by fake_advisor_id year: gen prop_christ=sumchris1out/newtotal

by fake_advisor_id year: egen sum_sec=sum(secular_school)
by fake_advisor_id year: gen sumsec1out=sum_sec-secular_school
by fake_advisor_id year: gen prop_sec=sumsec1out/newtotal

*Generate Interaction Terms accross religious schools
gen christianchrist2=christ_school*prop_christ
gen christianislam2=christ_school*prop_islam
gen christiansec2=christ_school*prop_sec

gen islamicchrist2=islam_school*prop_christ
gen islamicislam2=islam_school*prop_islam
gen islamicsec2=islam_school*prop_sec

gen secularchrist2=secular_school*prop_christ
gen secularislam2=secular_school*prop_islam
gen secularsec2=secular_school*prop_sec



*Standardize treatment variables*
egen stmusl=std(prop_islam)
egen stcath=std(prop_christ)
egen stsec=std(prop_sec)
 
gen christianchrist3=christ_school*stcath
gen  islamicchrist3=islam_school*stcath
 
gen christianislam3=christ_school*stmusl
gen islamicislam3=islam_school*stmusl
 
 *Gender of student
  tab StudGender, gen(sex)
 *Gender of advisor
  tab ADV_GENDER, gen(advsex)
 *Rank of advisor
  replace Rank="NA" if Rank==""  
  tab Rank, gen(advrank)
  
  *fill in missing SAT*
  sort ID, stable
	by ID: egen satmat = mean(math_sat)
	by ID: egen satverb = mean(verbal_sat)
	
*Create quality of peers variable
sort year fake_advisor_id, stable
by year fake_advisor_id: egen qualchris=mean(satmat) if christ_school ==1
by year fake_advisor_id: egen qualitychris=mean(qualchris)
gen christ_qualchris = qualitychris*christ_school
gen musl_qualchris = qualitychris*islam_school

by year fake_advisor_id: egen qualmusl=mean(satmat) if islam_school==1
by year fake_advisor_id: egen qualitymusl=mean(qualmusl)
gen christ_qualmusl = qualitymusl*christ_school
gen musl_qualmusl = qualitymusl*islam_school


by year fake_advisor_id: egen qualsec=mean(satmat) if secular_school==1
by year fake_advisor_id: egen qualitysec=mean(qualsec)
gen christ_qualitysec = qualitysec*christ_school
gen musl_qualitysec = qualitysec*islam_school

*Create gender of peers variable
gen sex = 1 if sex1==1
replace sex= 0 if sex2==1
by year fake_advisor_id: egen meanpeersex = mean(sex)
gen christmeanpeersex = christ_school* meanpeersex 
gen muslmeanpeersex = islam_school*meanpeersex

keep if year>=2002 & year<=2017 & bla==1 & FR_nonFR_status =="never enrolled as FR" &(prop_islam!=0 &  prop_christ!=0) & satmat!=. & satverb!=. & SemGPA!=. & MAJR_CODE!="MJRL" & qualitychris!=. & qualitymusl!=. 

tempfile advisor_chars
preserve
keep Rank Dept ADV_GENDER year fake_advisor_id
duplicates drop
sort fake_advisor_id year, stable
duplicates tag fake_advisor_id year, gen(dup)
drop if dup==1 & Rank=="NA"
collapse (first) ADV_GENDER Rank Dept, by( fake_advisor_id year)
save `advisor_chars',replace
restore


matrix A=0,0,0,0,0,0,0,0
set seed 11111993
encode Dept, gen(department)
forvalues year = 2002(1)2017 {
   display `year'
   quietly levelsof department, local(departments)
   foreach y of local departments{
   display "`y'"
   preserve
   quietly drop if year != `year'
   quietly drop if department != `y'
   quietly drop if christ_school ==.
   quietly drop if islam_school ==.
   quietly drop if secular_school ==.
   quietly drop if math_sat ==.
   quietly drop if verbal_sat ==.
   capture{	
   gen select = 0
   gen udraw = 0
	quietly levelsof fake_advisor_id, local(levels)
	foreach x of local levels {
	   noisily display `x'
 	   quietly replace select = 0
	   quietly replace select = 1 if fake_advisor_id== `x'
	   quietly summarize christ_school if select
	   scalar nsec = r(N)
	   scalar sumorig = r(sum)
	   scalar accum = 0
	   quietly summarize islam_school if select
	   scalar nsec2 = r(N)
	   scalar sumorig2 = r(sum)
	   scalar accum2 = 0
	   quietly summarize secular_school if select
	   scalar nsec3 = r(N)
	   scalar sumorig3 = r(sum)
	   scalar accum3 = 0
	   quietly summarize math_sat if select
	   scalar nsec4 = r(N)
	   scalar sumorig4 = r(sum)
	   scalar accum4 = 0	  
	   quietly summarize verbal_sat if select
	   scalar nsec5 = r(N)
	   scalar sumorig5 = r(sum)
	   scalar accum5 = 0
	   forvalues i = 1(1)10000 {
		quietly replace udraw = uniform()
		sort udraw, stable
		quietly summarize christ_school if _n <= nsec
		scalar sumfake=r(sum)
		scalar accum = accum + (sumfake < sumorig)  
		quietly summarize islam_school if _n <= nsec2
		scalar sumfake2=r(sum)
		scalar accum2 = accum2 + (sumfake2 < sumorig2) 
		quietly summarize secular_school if _n <= nsec3
		scalar sumfake3=r(sum)
		scalar accum3 = accum3 + (sumfake3 < sumorig3) 
		quietly summarize math_sat if _n <= nsec4
		scalar sumfake4=r(sum)
		scalar accum4 = accum4 + (sumfake4 < sumorig4)
		quietly summarize verbal_sat if _n <= nsec5
		scalar sumfake5=r(sum)
		scalar accum5 = accum5 + (sumfake5 < sumorig5) 
		}
		
		matrix A=A\ `year', `y' , `x',(accum/10000),(accum2/10000),(accum3/10000),(accum4/10000),(accum5/10000)
	}
	}
	restore
	}
}

svmat A
keep A1 A2 A3 A4 A5 A6 A7 A8
drop if _n==1
rename A1 year
rename A2 department
rename A3 advisorid
rename A4 pvalue
rename A5 pvalue2
rename A6 pvalue3
rename A7 pvalue4
rename A8 pvalue5
drop if year>2017



capture program drop KsChiSq
program KsChiSq 
 	matrix result=.,.,.
 forvalues year = 2002(1)2017 {
    quietly levelsof department, local(departments)
   foreach y of local departments{
     preserve
	quietly drop if year != `year'
	   quietly drop if department != `y'
	quietly drop if `1' == 0 | `1' == 1  	// To eliminate infeasible tests 
		  scalar teststat = 0
		  quietly sum `1'
		  scalar nobs = r(N)
		  forvalues i = 0.1(0.1)1 {
		    quietly sum `1' if ((`1' >= (`i'-0.1)) & (`1' < `i'))
		    scalar teststat = teststat + ((r(N)-0.1*nobs)^2)/(nobs*0.1)
		   }
		local pval = chi2(9,teststat)
		  display "`year', `y', `pval'"
		  matrix result=result\ `year', `y', `pval'
		restore
		}
 }
end

KsChiSq pvalue
svmat result
rename result1 res1 
rename result2 res2
rename result3 testpval1

sum testpval1
scalar christ_tests=r(N)
count if testpval1<0.05
scalar christ_fail=r(N)

KsChiSq pvalue2
svmat result
rename result3 testpval2
drop result*

sum testpval2
scalar islam_tests=r(N)
count if testpval2<0.05
scalar islam_fail=r(N)


KsChiSq pvalue3
svmat result
rename result3 testpval3
drop result*

sum testpval3
scalar sec_tests=r(N)
count if testpval3<0.05
scalar sec_fail=r(N)


KsChiSq pvalue4
svmat result
rename result3 testpval4
drop result*

sum testpval4
scalar math_tests=r(N)
count if testpval4<0.05
scalar math_fail=r(N)



KsChiSq pvalue5
svmat result
rename result3 testpval5
drop result*

sum testpval5
scalar verb_tests=r(N)
count if testpval5<0.05
scalar verb_fail=r(N)


rename advisorid fake_advisor_id
merge 1:1 fake_advisor_id year using `advisor_chars'
keep if _merge==3
gen female=( ADV_GENDER=="F")
gen tenure=(Rank=="Associate Professor" | Rank=="Chairperson & Associate Pr" | Rank=="Chairperson & Professor" | Rank=="Professor")
reghdfe pvalue female tenure, absorb(year Dept) cluster( fake_advisor_id)
scalar christ_fem=_b[female]
scalar christ_fem_se=_se[female]
scalar christ_ten=_b[tenure]
scalar christ_ten_se=_se[tenure]
reghdfe pvalue2 female tenure, absorb(year Dept) cluster( fake_advisor_id)
scalar islam_fem=_b[female]
scalar islam_fem_se=_se[female]
scalar islam_ten=_b[tenure]
scalar islam_ten_se=_se[tenure]
reghdfe pvalue3 female tenure, absorb(year Dept) cluster( fake_advisor_id)
scalar sec_fem=_b[female]
scalar sec_fem_se=_se[female]
scalar sec_ten=_b[tenure]
scalar sec_ten_se=_se[tenure]
reghdfe pvalue4 female tenure, absorb(year Dept) cluster( fake_advisor_id)
scalar math_fem=_b[female]
scalar math_fem_se=_se[female]
scalar math_ten=_b[tenure]
scalar math_ten_se=_se[tenure]
scalar obs=e(N)
reghdfe pvalue5 female tenure, absorb(year Dept) cluster( fake_advisor_id)
scalar verb_fem=_b[female]
scalar verb_fem_se=_se[female]
scalar verb_ten=_b[tenure]
scalar verb_ten_se=_se[tenure]




capture file close TableA2
file open TableA2 using "$dir\Tables\TableA2.tex",write replace
file write 	TableA2 /// 
"\begin{tabular}{l*{6}{c}}" _n ///	
	"\toprule" _n ///	
	"&\multicolumn{1}{c}{(1)}&\multicolumn{1}{c}{(2)}&\multicolumn{1}{c}{(3)}&\multicolumn{1}{c}{(4)}&\multicolumn{1}{c}{(5)}\\" _n ///	
	"&\multicolumn{1}{c}{  Math SAT  }&\multicolumn{1}{c}{  Verbal SAT  }&\multicolumn{1}{c}{  secular School  } &\multicolumn{1}{c}{  Christian School   }&\multicolumn{1}{c}{  Islam School  }\\ " _n ///	
	"&\multicolumn{1}{c}{  Empirical   }&\multicolumn{1}{c}{ Empirical   }&\multicolumn{1}{c}{  Empirical   }&\multicolumn{1}{c}{ Empirical }&\multicolumn{1}{c}{  Empirical }\\" _n ///	 
	"&\multicolumn{1}{c}{  P-Value   }&\multicolumn{1}{c}{ P-Value   }&\multicolumn{1}{c}{  P-Value   }&\multicolumn{1}{c}{ P-Value }&\multicolumn{1}{c}{  P-Value }\\" _n ///	
	 "\midrule" _n ///	
	 "\addlinespace" _n ///	
	 "\multicolumn{1}{c}{\textbf{A. Test for Student Characteristics  }} \\" _n ///	
	 "\addlinespace" _n ///	
	 "\addlinespace" _n ///	
	 "$\chi^2$ goodness of fit test          & "      %8.0gc (math_fail) "/"  %8.0gc (math_tests) "& "    %8.0gc  (verb_fail) "/"  %8.0gc (verb_tests) "& "  %8.0gc (sec_fail) "/"  %8.0gc (sec_tests) "& "  %8.0gc (christ_fail) "/"  %8.0gc (christ_tests) "& "  %8.0gc (islam_fail) "/"  %8.0gc (islam_tests)  "\\" _n ///	 
	 "(no. failed/total tests)  \\" _n ///	
	 "\addlinespace" _n ///	
	 "\midrule" _n ///	
	 "\addlinespace" _n ///	
	 "\multicolumn{1}{c}{\textbf{B. Test for Advisor Characteristics  }} \\" _n ///	
	 "\addlinespace" _n ///	
	 "Female Advisor &" %9.3f (math_fem)  "&" %9.3f (verb_fem)  "&" %9.3f (sec_fem)  "&" %9.3f (christ_fem)  "&" %9.3f (islam_fem)  "\\" _n ///	
	 "&" "(" %-9.3f (math_fem_se) ")" "&" "(" %-9.3f (verb_fem_se) ")" "&" "(" %-9.3f (sec_fem_se) ")" "&" "(" %-9.3f (christ_fem_se) ")" "&" "("%-9.3f (islam_fem_se) ")" "\\" _n ///	
	 "\addlinespace" _n ///	
	 "Associate/Full Professor &" %9.3f (math_ten)  "&" %9.3f (verb_ten)  "&" %9.3f (sec_ten)  "&" %9.3f (christ_ten)  "&" %9.3f (islam_ten)  "\\" _n ///	
	 "&" "(" %-9.3f (math_ten_se) ")" "&" "(" %-9.3f (verb_ten_se) ")" "&" "(" %-9.3f (sec_ten_se) ")" "&" "(" %-9.3f (christ_ten_se) ")" "&" "(" %-9.3f (islam_ten_se) ")" "\\" _n ///
	 "\(N\)       &"     %8.0gc (obs)      "& " %8.0gc (obs)   "& " %8.0gc (obs)   "& " %8.0gc (obs)   "& " %8.0gc (obs)   "\\" _n ///
	 "\addlinespace" _n ///
	 "\toprule" _n ///			
		"\end{tabular}\\"								
file close TableA2	
